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Abstract 

00 \ The initial gluon multiplicity per unit area per unit rapidity, dN jl? jdri, in high 

O ■ energy nuclear collisions, is equal to j ^{^g^ ixL^g^ iif' j g^ , with \i^ proportional to 



l> 



the gluon density per unit area of the colliding nuclei. For an SU(2) gauge theory, 
^^ \ we compute f]ss{g^[xV) = 0.14 it 0.01 for a wide range in g'^fiL. Extrapolating 

O I to SU(3), we predict dN/L"^ /drj for values of g^iJ.L in the range relevant to the 

-««^ ■ Relativistic Heavy Ion Collider and the Large Hadron Collider. We compute the 

'^ initial gluon transverse momentum distribution, dN/L'^/d'^k±, and show it to be 






X 



well behaved at low k±. 

A topic of considerable current interest is the possibility of forming an equilibrated 
plasma of quarks and gluons, a quark-gluon plasma (QGP), in very high energy nuclear 
collisions. Experimental signatures of such a plasma may provide insight into the nature 
5^ I of the QCD phase diagram at finite temperature and baryon density ||l|. 

Equally interesting, is the information that heavy ion collisions may provide about 
the distributions of partons in the wavefunctions of the nuclei before the collision. At 
very high energies, the growth of parton distributions in the nuclear wavefunction satu- 
rates, forming a state of matter sometimes called a Color Glass Condensate (CGC) 0. 
The condensate is characterized by a bulk momentum scale Qg. If Qs ^ ^qcd, the 
properties of this condensate, albeit non-perturbative, can be studied in weak coupling. 
The partons that comprise this condensate are freed in a collision. Since the scale of 
the condensate, Qg, is the only scale in the problem, the initial multiplicity and energy 
distributions of produced gluons at central rapidities are determined by this scale alone. 
We will briefly discuss later the relation of these quantitities to physical observables. 

The above statements may be quantified in a classical effective field theory approach 
(EFT) to high energy scattering [Q. The EFT is classical because, at central rapidities, 
where a; <^ 1 and p± ^ ^qcd, {x ~ P±/a/s), parton distributions grow rapidly with 



decreasing x giving rise to large occupation numbers. Briefly, the EFT separates partons 
in a liadron (or nucleus) into static, high x valence and hard glue sources, and "wee" 
small X fields. For a large nucleus in the infinite momentum frame, the hard sources 
with color charge density p, are randomly distributed in the transverse plane with the 
distribution 

Pi[p]) = exp (-^ Jd'x^p\x^)^ . (1) 

The average squared color charge per unit area is determined by the parameter /x^, 
which is the only dimensional parameter in the EFT apart from the linear size L of the 
nucleus. Parton distributions, correlation functions of the wee gauge fields, are computed 
by averaging over gauge fields with the weight P{[p]) [§. 

Quantum corrections to the EFT P are implemented using Wilson renormalization 
group techniques [0. The scale /x^ grows with decreasing x, and can be estimated from 
the nucleon quark and gluon distributions at x <^ 1 |^. The saturation scale — at which 
parton distributions stop growing rapidly with decreasing x — is Q^ ~ Gg'^fi/An 0, a 
function determined self-consistently from the typical x and Q^ of interest. Since most of 
the saturated partons have momenta of this order, Qs, rather than g'^fi, is the appropriate 
scale. At RHIC energies, Q^ ~ 1 GeV, and at LHC, Qs ~ 2-3 GeV. The magnitude of 
Qs is a rough estimate because of our uncertain knowledge of gluon structure functions 
at these energies. This scale, and the properties of the CGC, may be further determined 
in future high energy deeply inelastic scattering (DIS) experiments off nuclei [^. 

In this letter, we obtain a non-perturbative relation between the multiplicity of glu- 
ons produced in heavy ion collisions at central rapidities on the one hand, and g'^fi, or 
equivalently Qs, on the other. We will also demonstrate that non-perturbative strong 
field effects, at momenta kj_ ^ Qs, qualitatitively alter transverse momentum distribu- 
tions, rendering them infrared finite. A preliminary version of these results was reported 
in Ref. [Q. In a previous letter [^, we obtained a similar expression for the energy of 



gluons, per unit rapidity, produced shortly after a very high energy nuclear collision. 
The problem of initial conditions [|T6| for nuclear scattering can be formulated in the 



classical EFT |]T0[ in the gauge A"^ = 0. Matching the Yang-Mills equations D^F'^'^ = J" 
in the four light cone regions, along the light cone, one obtains for the gauge fields 
in the forward light cone, at proper time r = 0, the relations A^ = A\ + A2 ^^^ 
A^ = ±.igx'^[A\, A^^/2. Here J^ = Si,25'''^5(2;^)p^'^(a^±) are random light cone sources 
corresponding to the valence or hard glue sources in the two nuclei. The transverse 
pure gauge fields A\2{p'^), with i = 1,2 are solutions of the Yang-Mills equations for 
each of the two nuclei before the collision. With these initial conditions, the Yang- 
Mills equations can be solved in the forward light cone to obtain gluon configurations 
at late proper times. Since the initial conditions depend on the sources p^, averages 
over different realizations of the sources — specified by the weight in Eq. (|l|) — must be 
performed. 

Perturbative solutions for the number distributions in transverse momentum, per 
unit rapidity, were obtained in Refs. |]TD|, |T^. These were shown to be infrared divergent. 



In the classical EFT, this divergence is logarithmic. The number distributions have the 
form 

as \ k± J \asiij 

for k±_ S> asiJi- The perturbative description breaks down at /c^ ~ 05/1. Thus, for 
robust predictions of gluon multiplicity distributions, a fully non-perturbative study of 
the classical EFT is necessary ||12|| . 

The model is discretized on a lattice in the transverse momentum plane. Boost in- 
variance and periodic boundary conditions are assumed. The lattice Hamiltonian is the 
Kogut-Susskind Hamilitonian in 2+1-dimensions coupled to an adjoint scalar field. The 
lattice field equations are then solved by computing the Poisson brackets, with initial 
conditions that are the lattice analogs of the continuum initial conditions mentioned ear- 
lier. Technical details of our simulations can be found in Refs. [0, |T^. Our simulations 



are presently only for an SU(2) gauge theory — the full SU(3) case will be studied later. 

The scale g'^^ and the linear size of the nucleus L are the only physically relevant 
dimensional parameters of the classical EFT. Any dimensional quantity P well defined 
within the EFT can then be written as (g'^fi)'^ fp{g'^fiL), where d is the dimension of 
P. All the non-trivial physical information is therefore contained in the dimensionless 
function fp{g'^fiL). On the lattice, P will generally depend on the lattice spacing a; this 
dependence can be removed by taking the continuum limit a — *^ 0. Assuming g = 2, 
the physically relevant values of g'^fi for RHIC and LHC energies are ~ 2 GeV and ~ 4 
GeV respectively. Also, assuming central Au-Au collisions, we obtain L = 11.6 fm as the 
physical linear dimension of our square lattice. Thus, g'^fJ^L ^ 120 (240) for RHIC (LHC). 
In Ref. [jl5[, we computed the initial energy density, per unit area per unit rapidity, to 
be dE / L? jdrj = /e {g'^jj)'^. We computed Je as a function of g^fJ^L, and extrapolated 
our results for /e to the continuum limit. We could however just as well expressed our 
result as dE/nR'^/dri = ce{QsR) {Qs)^^ with ce ~ 4.3-4.9 in the region of interest. 

We will now report on our results for the initial multiplicity of gluons produced at 
central rapidities in very high energy nuclear collisions. This quantity, while not directly 
observable, is related to the number of hadrons produced at central rapidities [^]. The 
initial multiplicity and momentum distribution of gluons also determine the equilibration 
time, the temperature and the chemical potential of the QGP [^, ^. The various 



signatures of QGP formation are highly sensitive to these quantitites |[^. 

We must first clarify what we mean by the number of quanta in the interacting, non- 
Abelian gauge theory at hand. To motivate the discussion, let us first consider a free 
field theory whose Hamiltonian in momentum space has the form 

Hf = \Y.{Hk)?+^\m<pm'). (3) 

k 

where 0(A;) is the fcth momentum component of the field, '^{k) is its conjugate momentum, 
and ui{k) is the corresponding eigenfrequency. The average particle number of the /c-th 



mode is then 

Nik) = ujik)mm = ^/mmmf), ■ (4) 

In our case, the average () is over the initial conditions. 

Clearly, any extension of this notion to interacting theories should reduce to the 
standard free- field definition of the particle number in the weak coupling limit. However, 
this requirement alone does not define the particle number uniquely outside a free theory. 
We therefore use two different generalizations of the particle number to an interacting 
theory. Each has the correct free-field limit. Even though the fields in question are 
strongly interacting at early times, they are only weakly coupled at late times, and it is 
only at this stage that it becomes reasonable to define particle number. We verify that 
the two definitions agree in this weak-coupling regime. 

Our first definition of the multiplicity is straightforward. We impose the Coulomb 
gauge condition in the transverse plane, V± ■ A±_ = 0, and substitute the momentum 
components of the resulting field configuration into Eq. (p. One option now is to as- 
sume uj{kj_) to be the standard massless (lattice) dispersion relation and use the middle 
expression of Eq. (^) to compute N{k^). Alternatively, we can determine N{k±) from 
the rightmost expression of Eq. (^; the middle expression of Eq. (^) can then be used 
to obtain u!{k±). The second option is preferable; it does not require us to assume that 
the dispersion relation is linear. 

Our second definition is based on the behavior of a free-field theory under cooling. 
Consider a simple relaxation equation for a field in real space, 

dt(f){x) = -dH/d(j){x), (5) 

where t is the cooling time [not to be confused with real or proper time) and H is the 
Hamiltonian. For a free field theory {H = Hf) the relaxation equation has exactly the 
same form in the momentum space with the solution 0(/c, t) = (j){k, 0) exp(— a;^(/c)t). The 
potential energy of the relaxed free field is V{t) = {l/2)J2k^'^{k)\(j){k,t)\'^. It is then 
easy to derive the following integral expression for the total particle number of a free-field 
system: 

Eq. d^) can be solved numerically for interacting fields. Subsequently, V{t) can be deter- 
mined, and A^ can be computed by numerical integration. Note that in a gauge theory 
the relaxation equations are gauge-covariant, and the relaxed potential V{t) is gauge- 
invariant, entailing gauge invariance of this definition of the particle number. This is an 
attractive feature of the cooling method. On the other hand, unlike the Coulomb gauge 
computation discussed earlier, this cooling technique presently only permits determina- 
tion of the total particle number. It cannot be used to find the number distribution 
N{k^). 

Both our definitions cease to make sense if the system is far from linearity. In partic- 
ular, if the theory has metastable states, and the system relaxes to one of these states, 



the right hand side of (|^) will diverge. The expression in Eq. (^ can then still be for- 
mally determined in the Coulomb gauge but its interpretation as a particle number is 
problematic. This situation deserves special attention and will be discussed in detail 
elsewhere. We did not observe any effects of metastabilty within the range of parame- 
ters of the current numerical study. In particular, we verified the convergence of Eq. (^ 
with respect to the upper limit of integration. 

We now present our results using both the techniques discussed. We begin with the 
number distribution, which can only be computed in Coulomb gauge. We have verified, 
for g^iiL = 35.35, that in the range of values of g'^fia considered here the system is close 
to the continuum limit. This is consistent with our earlier analysis of the lattice spacing 
dependence of a more ultraviolet-sensitive quantity, the energy density |T^. In Fig. la, 
we plot the gluon number distribution, n{k±) = dN/L? /dk^ = iV(fcj_)/(27r)^ versus k± 
for fixed g^fJ^L = 35.5, but for different values of the lattice spacing g'^fia. For large 
kj_, the finest lattice {g'^fJ^a = 0.138) agrees well with the lattice perturbation theory 
(LPTh) analogue of Eq. (|^). At smaller k±, the distribution is softer, and converges to a 
constant value. In Fig. lb, we plot the gluon distribution in the infrared, at fixed g'^fia, 
for different g^fiL (148.5 and 297). We notice that these distributions are nearly universal 
and independent of g^^Ll Also, the convergence of the distribution to a constant value 
is more clearly visible in Fig. lb. 

When k± < g^fi, non-perturbative effects qualitatively alter the perturbative num- 
ber distribution, rendering it finite in the infrared. Unfortunately, since these effects 
are large, an analytical understanding of the behavior at low k± is lacking. Our re- 
sults, despite being universal, are not simply fit by any of the physically motivated 
parametrizations we have considered. 

From our previous discussion, the formula 

^f-?^"*'^''^)''^. (7) 

relates the number of produced gluons per unit area per unit rapidity at zero rapidity to 
g'^fi. We have computed /at on the lattice using the two techniques discussed earlier. Our 
result for f^ as a function of g^fiL, for the smallest values of g^fia feasible, are plotted 
in Fig. 2. We see that the agreement between the cooling and Coulomb techniques at 
larger values of g^fiL is excellent. It is not as good at the smaller values - in general, 
the cooling number is more reliable [113. We also note that Fig. 2 demonstrates that the 



distributions in Fig. lb are not quite universal — otherwise, /^v would be a constant. 
We see instead that it has a weak logarithmic rise with g^fiL for larger g'^iiUs. Table 1 
lists /at for various g^jiL. The third row is the Coulomb gauge number after cooling — 
seeRef. [|l|]. 

In Fig. 3, we plot the dispersion relation uj{k±) vs k\_ using the relation Eq. (^. 
All the dispersion curves rapidly approach the uj{k±) = k± asymptote characteristic of 
on-shell partons, while exibiting a mass gap at zero momentum. We reserve for a later 
work a detailed study of this mass gap and its role in rendering the number distributions 
infrared finite. 



g^^L 


35.36 


70.71 


106.1 


148.5 


212.1 


297.0 


g^^ia 


.276 


.276 


0.207 


.29 


.41 


.29 


JN (cooling) 


.116 ±.001 


.119 ±.001 


.127 ±.001 


.138 ±.001 


.146 ±.001 


.151 ±.001 


/at (Coulomb) 


.127 ±.002 


.125 ± .002 


.135 ±0.001 


142 ± .001 


.145 ±.001 


.153 ±.001 


fN (res.) X 103 


14 ±2 


7.8 ±0.2 


8.9 ±0.2 


5.6 ±0.1 


7.12 ±0.08 


4.83 ± .04 



Table 1: Values of f^ vs g'^fiL, for fixed g^fia, plotted in Fig. 0. fjy (res.) is defined in 
Ref. ^ 



We can compare our results for the number distribution to the one predicted by A. 
H. Mueller [^]. In terms of Qg and R, we can re-write Eq. |^ as 

1 dN N^-l 1 ^2 

c^ -^T — T-9 — Qs ■ 



71 R^ drj Nc An'^as 

Mueller estimates the non-perturbative coefficient c^ to be of order unity. If we take 
/tv = 0.14 ± 0.01, as is the case for much of the range studied, we find cat = 1.29 ± 0.09, 
a number of order unity as predicted by Mueller. Despite this close agreement, the 
transverse momentum distributions, shown in Fig. la and lb, and discussed earlier, look 
quite different from Mueller's guess of 0{Ql — kj_). The ^-function distribution was only 
a rough guess to represent a qualitative change of the distributions at k± ^^ Qg. 

A large number of models of particle production in nuclear collisions at RHIC and 
LHC energies can be found in the literature. A nice recent summary of their various 
predictions and relevant references can be found in the compilation of Ref. [^ . Naively 
extrapolating our results to SU(3), we find for Au-Au central collisions at RHIC energies, 
dN/dr] ~ 950 for /^ = 0.132 ± 0.006 {g'^^iL ^ 120 — we take the mean of the 106 and 
148 cooling point). Similarly, for LHC energies f^ = 0.148 ± 0.002 {g'^fiL ~ 255 — 
the mean of the 212 and 297 cooling point), one finds dN/drj ~ 4300. In particular, 
comparing our predictions with those of pQCD based models |]19|], we find our numbers 
to be in rough agreement. However, if we include a K factor like many of these models 
do, our numbers will be roughly a factor of 2 larger. 

There is considerable uncertainity in the value of Qs because the gluon densities at 
the relevant x and Q^ are ill-known. Since the multiplicity depends quadratically on 
Qs, a prediction of the same is perforce unreliable. Distinguishing between different 
models will therefore require, at the very least, testing their predictions for the scaling 
of multiplicities with A and with ^/s. In our case, Qs ~ A^^^, hence from Eq. |^ the 
number per unit rapidity will, up to logarithms of A, be proportional to A. In the 
EFT, one naively expects the energy dependence to be a power law, the power being 
determined by the rise in the gluon density at small x. Quantitative estimates of the 
dependence of Qs with energy in the saturation region are being developed. Predictions 
from other models vary significantly, ranging from a power law dependence [^ of Qs to 



Qs oc exp(-i/ln(s/so)), where Sq is a constant. In the latter case [23|, it is claimed that a 
good fit to the multiplicity from existing high energy hadron scattering data is obtained. 
Data from RHIC will help constrain the energy dependence of Qs- 

6 



The reader should also note that our relations are derived only for the initial parton 
multiphcity distributions at central rapidities. These provide the initial conditions for 
the subsequent evolution of the system, which can be investigated in a transport ap- 
proach [jl^ ^ . The rate of chemical equilibration, and uncertainities due to hadroniza- 
tion also have to be taken into account in predictions of observables such as charged 
hadron multiplicities. Conversely, these observables may help constrain the saturation 
scale Qsi and inform us about the very earliest stages of nuclear collisons. 

In summary, we have derived a non-perturbative relation between the multiplicity 
of produced partons and the saturation scale of parton distributions in high energy 
nuclear collisions. We have computed number distributions which have the predicted 
perturbative behavior in the ultraviolet, and are finite in the infrared. At present, in our 
approach, we are only able to make qualitative "ball-park" predictions. However, we have 
developed a framework in which these can be quantified and extended in a consistent 
manner to study a large number of final state observables in heavy ion collisions. 
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Figure 1: a: n{k^) = dN/ L?' /(Pk± as a function of the gluon momentum k for g^^L = 
35.35 and the values 0.138 (squares), 0.276 (plusses), and 0.552 (diamonds) oi g^fta. The 
gluon momentum k is in units of g'^fi. The solid line is a fit of the lattice analog of the 
perturbative expression Eq. (|^) to the high- momentum part of the g^fia = 0.138 data, b: 
n{k±) at soft momenta at g^fia = 0.29 for the values 148.5 (plusses) and 297 (diamonds) 
of g'^i^L. 
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Figure 2: The function /^r, defined in Eq. (|^ as a function of g^fiL, obtained by the 
relaxation method (plusses) and by the Coulomb gauge fixing (diamonds). The values of 
g'^fxa are 0.276 for g^fiL = 35.35 and g^fxL = 70.8; 0.29 for g^fiL = 148.5 and g^fiL = 297; 
and 0.414 for g^fxL = 212. 




Figure 3: Gluon dispersion relation uj{k±) obtained from Eq. @, for the values 70. 
(diamonds), 148.5 (plusses), and 297 (squares), with the values of g^fia as in Figure 0. 
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